# Gilardi, Fabrizio, "The Temporary Importance of Role Models for Women's Political Representation", American Journal of Political Science
# Code to replicate Figure SI6 (Effect of being ten minutes closer to the next women’s group on the number of female candidates)
# gilardi@ipz.uzh.ch, 2014-06-24

# Set working directory
setwd("../Data/")

# Load packages
library(Zelig)

# Load data
d <- read.csv("dataset-full.csv")

# Function to compute mode
Mode <- function(x) {
 ux <- unique(x)
  ux[which.max(tabulate(match(x, ux)))]
}


# Estimate models

out.2.z <- glm(n.women.cand ~ sl.autob.n.women.elected + n.women.elected.lagged + n.women.cand.previously.elected + n.men.cand.previously.elected + n.women.cand.lagged + n.people.elected + n.cand + referenda + wom.kr + log.dist.zh + log.pop + kr.svp + steuerfuss + agglo + land + dist.fp
	+ dist.fp*y1978
	+ dist.fp*y1982
	+ dist.fp*y1986
	+ dist.fp*y1990
	+ dist.fp*y1994
	+ dist.fp*y1998
	+ dist.fp*y2002
	+ dist.fp*y2006
	+ dist.fp*y2010
	, data=d, family="poisson")
summary(out.2.z)

mm <- model.frame(out.2.z)

z.out.2 <- zelig(n.women.cand ~ sl.autob.n.women.elected + n.women.elected.lagged + n.women.cand.previously.elected + n.men.cand.previously.elected + n.women.cand.lagged + n.people.elected + n.cand + referenda + wom.kr + log.dist.zh + log.pop + kr.svp + steuerfuss + agglo + land + dist.fp
	+ dist.fp*y1978
	+ dist.fp*y1982
	+ dist.fp*y1986
	+ dist.fp*y1990
	+ dist.fp*y1994
	+ dist.fp*y1998
	+ dist.fp*y2002
	+ dist.fp*y2006
	+ dist.fp*y2010
	, model="poisson", data=mm)
summary(z.out.2)


# Vector of election years
yr <- seq(1974, 2010, 4)

# Significance level
sig <- qnorm(0.975)

# Empty matrix to store simulated values
pr.l.2 <- pr.h.2 <- pr.d.2 <- matrix(NA, length(yr), 4)
pr.l.2[,1] <- pr.h.2[,1] <- pr.d.2[,1] <- yr


# Simulate values and store them in matrix

for(i in 1:length(yr)){

# Baseline distance from women's group in minutes
fp.l <- 20
# Reduce distance by ten minutes
fp.h <- fp.l - 10

# Set values for longer distance
x.l <- setx(z.out.2, sl.autob.n.women.elected=round(Mode(d$sl.autob.n.women.elected[d$year==yr[i]]), 2), n.women.elected.lagged=round(Mode(d$n.women.elected.lagged[d$year==yr[i]]), 1), n.women.cand.lagged=round(Mode(d$n.women.cand.lagged[d$year==yr[i]]), 1), n.women.cand.previously.elected=round(Mode(d$n.women.cand.previously.elected[d$year==yr[i]]), 1), n.men.cand.previously.elected=round(Mode(d$n.men.cand.previously.elected[d$year==yr[i]]), 1), dist.fp=fp.l,
	y1978 = unique(d$y1978[d$year==yr[i]]),
	y1982 = unique(d$y1982[d$year==yr[i]]),
	y1986 = unique(d$y1986[d$year==yr[i]]),
	y1990 = unique(d$y1990[d$year==yr[i]]),
	y1994 = unique(d$y1994[d$year==yr[i]]),
	y1998 = unique(d$y1998[d$year==yr[i]]),
	y2002 = unique(d$y2002[d$year==yr[i]]),
	y2006 = unique(d$y2006[d$year==yr[i]]),
	y2010 = unique(d$y2010[d$year==yr[i]])
	)

# Set values for shorter distance
x.h <- setx(z.out.2, sl.autob.n.women.elected=round(Mode(d$sl.autob.n.women.elected[d$year==yr[i]]), 2), n.women.elected.lagged=round(Mode(d$n.women.elected.lagged[d$year==yr[i]]), 1), n.women.cand.lagged=round(Mode(d$n.women.cand.lagged[d$year==yr[i]]), 1), n.women.cand.previously.elected=round(Mode(d$n.women.cand.previously.elected[d$year==yr[i]]), 1), n.men.cand.previously.elected=round(Mode(d$n.men.cand.previously.elected[d$year==yr[i]]), 1), dist.fp=fp.h,
	y1978 = unique(d$y1978[d$year==yr[i]]),
	y1982 = unique(d$y1982[d$year==yr[i]]),
	y1986 = unique(d$y1986[d$year==yr[i]]),
	y1990 = unique(d$y1990[d$year==yr[i]]),
	y1994 = unique(d$y1994[d$year==yr[i]]),
	y1998 = unique(d$y1998[d$year==yr[i]]),
	y2002 = unique(d$y2002[d$year==yr[i]]),
	y2006 = unique(d$y2006[d$year==yr[i]]),
	y2010 = unique(d$y2010[d$year==yr[i]])
	)


# Simulate values for longer distance
s.out.l <- sim(z.out.2, x=x.l)
# Simulate values for shorter distance
s.out.h <- sim(z.out.2, x=x.h)
# Simulate first differences
s.out.d <- sim(z.out.2, x=x.l, x1=x.h)

# Expected value for longer distance
pr.l.2[i,2] <- summary(s.out.l)$stats$"Expected Values: E(Y|X)"[1]
# Lower bound of confidence interval for longer distance
pr.l.2[i,3] <- pr.l.2[i,2] - sig*summary(s.out.l)$stats$"Expected Values: E(Y|X)"[2]
# Upper bound of confidence interval for longer distance
pr.l.2[i,4] <- pr.l.2[i,2] + sig*summary(s.out.l)$stats$"Expected Values: E(Y|X)"[2]


# Expected value for shorter distance
pr.h.2[i,2] <- summary(s.out.h)$stats$"Expected Values: E(Y|X)"[1]
# Lower bound of confidence interval for shorter distance
pr.h.2[i,3] <- pr.h.2[i,2] - sig*summary(s.out.h)$stats$"Expected Values: E(Y|X)"[2]
# Upper bound of confidence interval for shorter distance
pr.h.2[i,4] <- pr.h.2[i,2] + sig*summary(s.out.h)$stats$"Expected Values: E(Y|X)"[2]

# Expected value of first differences
pr.d.2[i,2] <- summary(s.out.d)$stats$"First Differences: E(Y|X1) - E(Y|X)"[1]
# Lower bound of confidence interval for first differences
pr.d.2[i,3] <- pr.d.2[i,2] - sig*summary(s.out.d)$stats$"First Differences: E(Y|X1) - E(Y|X)"[2]
# Upper bound of confidence interval for first differences
pr.d.2[i,4] <- pr.d.2[i,2] + sig*summary(s.out.d)$stats$"First Differences: E(Y|X1) - E(Y|X)"[2]

}


#pdf(file="Figure-SI6.pdf", paper="special", width=7, height=5.5)
par(mar=c(3.5,3.5,2,1), mgp=c(2.5,0.8,0), cex.axis=1.2, cex.lab=1.2, font.main=1)
plot(yr, pr.d.2[,2], type="n", ylim=range(pr.d.2[,2:4]), xlab="", ylab="Expected increase in nr of women cand.", axes=F)
axis(1, at=yr)
axis(2)
abline(h=0, lty=2)
points(yr, pr.d.2[,2], type="b", cex=1.5)
segments(yr, pr.d.2[,3], yr, pr.d.2[,4])
legend("topright", pch=c(), legend="Effect of being 10 minutes closer\nto the nearest women's group", pt.cex=1.5, cex=1.2, bty="n")
#dev.off()





